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ABSTRACT 

The light-curves of the two large amplitude variable stars, R Scuti and AC 
Herculis, display irregular pulsation cycles with 'periods' of 75 and 35 days, re- 
spectively. We review the evidence that the observed time-series are generated 
by low dimensional chaotic dynamics. In particular, a global flow reconstruc- 
tion technique indicates that the minimum embedding dimensions are 4 and 
3, for R Set and AC Her, respectively, whereas the fractal dimensions are in- 
ferred to be di 3.1 and 2.05^^^2.45, respectively. It thus appears that 
the dimensions of the dynamics themselves are also 4 and 3. 

1. Introduction 

For several centuries it has been known that a sizeable subset of the stars ex- 
hibits a variable energy output (luminosity). The best known and best studied of 
these variable stars are the Cepheids that have a clock-like periodicity. Their noto- 
riety comes from the tight relation between their period and luminosity and their 
consequent role as cosmological distance indicators. In parallel, there exists a class 
of variable stars with light-curves that instead exhibit various degrees of irregularity. 
These objects undergo changes in luminosity that can be as large as a factor of forty, 
indicating violent radial (i.e. t = 0) pulsations. Very little theoretical attention had 
been devoted to these stars. 

Why do we think that the dynamics underlying these irregular pulsations could 
be low dimensional and chaotic? One of the reasons is that decade old numerical 
hydrodynamical simulations showed low dimensional chaotic behavior and period 
doubling cascades in models.EH However, the observational confirmation of low di- 
mensional chaos in these stars, that is the subject of this paper, had to wait a little 
longer for a number of reasons. 

2. The observational data 

Nonlinear techniques are new in variable star Astronomy, and it is therefore 
hard to find data sets that are suitable for a modern nonlinear time series analy- 
sis and the search for low dimensional chaos. Partially because of the long 'periods' 



involved, typically tens to hundreds of days, rarely are enough pulsation cycles cov- 
ered, or they are covered incompletely. Generally, professional astronomers have not 
seen any merit in dedicating the necessary telescope time to what the consider the 
gathering of 'useless' erratic data. In contrast, thanks to the observations of numerous 
amateur astronomers, the American Association of Variable Star Observers (AAVSO) 
has recently made available long and relatively well sampled light-curve data of two 
stars, R Scuti and AC Herculis. (In Astronomy variable stars were named after the 
constellation in which they appear preceded by one or two capital letters, until the 
letters were exhausted). Both objects are characterized by irregular cycles of ~ 75 d 
and 35 d, respectively, with alternating shallow and deep minima. The details of the 
analyses of these light-curves have been published elsewherel™ Here we present a 
summary review. 
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Fig. 1 Typical segments of observational lightcurve data (dots) and fits (lines), magnitude versus 
time [years]; R Scuti (top) and AC Herculis (bottom). 

Fig. 1 displays typical subsections of the two data sets (points) with the solid 
line showing our smoothed and interpolated fit (see below). These data have been 
taken visually by many amateur observers worldwide at essentially random times. 
Therefore there is considerable scatter in the data and there are some observational 
gaps. We have plotted the lightcurves in magnitudes (for historical reasons defined as 
—2.5 Logx luminosity) rather than the more 'physical' luminosity. For both data sets 
the error in the magnitudes is found to have a normal distribution, independent of 
the magnitude (which is not true for the luminosity itself). The reason is to be found 
in the visual nature of the observations and the logarithmic response of the human 
eye. We thus work with the magnitude. Although the individual observations have a 
large error associated with them (0.20 and 0.15 mag, in R Set and AC Her, resp.) it is 
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possible to extract a reasonably accurate light-curve by averaging and interpolating. 
The individual errors are larger for R Set, but there are about four times as many 
data points and the pulsation amplitude is about a factor of 4 larger, so that the error 
on the final interpolated curve is much smaller for R Set than for AC Her. 
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Fig. 2 Observational lightcurves of R Scuti and best synthetic signal (4D). 
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Fig. 3 Observational lightcurves of AC Herculis (top) and best synthetic signal (5D). 

The preparation of the data is a delicate process. We want to eliminate as 
much observational and other extrinsic noise as possible. However there is additional 
'intrinsic noise' associated with turbulence and convection in the star that, strictly 
speaking, is part of the dynamics. It is well known that well developed turbulent 
motions constitute high-dimensional behavior, and a priori it may be objected this 
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would preclude us from uncovering any low dimensional dynamics. The success of 
our nonlinear analysis provides an a posteriori indication that the pulsation can be 
separated into a large amplitude 'true' pulsation, and a low amplitude pulsational 
jitter. Our pre-smoothing eliminates this high- dimensional jitter together with the 
true noise, leaving a signal that represents the 'true' low dimensional dynamics. 

To be specific, the bare data were first averaged in 5.0 day bins for R Set and in 
2.5 day bins for AC Her, followed by a cubic spline smoothing, where the smoothing 
parameter a was a free parameter in the flow reconstruction below. The result is then 
a lightcurve data set that is sampled at equal time-intervals. 

Figs. 2 and 3 in the top rows show the smoothed and interpolated observational 
data of R Set and AC Her, respectively. In Fig. 4 we display the (normalized) am- 
plitude Fourier spectra of the two observed lightcurves. The FS of R Set has much 
broader features than the more regular AC Her as is to be expected from the appear- 
ance of the lightcurves. 
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Fig. 4 Amplitude Fourier spectrum of R Set and AC Her lightcurve data. 



It is also of interest to look at the Broomhead-King (BK) projections^ of a 
signal. These projections onto the eigenvectors of the correlation matrix are optimal 
for spreading out the features of the signal. The lowest BK projections for R Set and 
AC Her are thus shown in columns 1 and 5 in Fig. 5. 

Before proceeding with a nonlinear analysis aimed at discovering and describing 
low dimensional chaos we briefly mention further tests to ascertain that the signals 
cannot be the result of a multi-periodic (n-torus) dynamics. 

For R Set we have performed the following tests J Taking the frequencies of the 
highest 35 peaks in the FS (Fig. 2) we have made a multi-periodic Fourier fit for 
the 35 amplitudes and phases. The fit is quite good over the temporal range of the 
data set, but when it is extrapolated in time the signal totally fails to resemble the 
R Set lightcurve. Similarly, a linear autoregressive process with white noise (AR or 
ARMA) bears little resemblance to the data. Both tests indicate that the process that 
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produces the lightcurve is strongly nonlinear .Q'Q 

For AC Her, the situation is not so clear-cut because of the high noise level. 
We note here prewhitening with the frequencies of the 24 highest FS peaks leaves a 
noise level that is at least 3-4 times the noise that is expected from the observational 
error.! 
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Fig. 5 Broomhead-King projections for the R Set and AC Her data, as well as for the respective best 
synthetic signals; from top to bottom, coordinates (£2,£i),(£3,£i)and (£3, £2)- 

There are theoretical reasons for not thinking that so many frequencies could 
be observed. First these stars are believed to be radial pulsators, and the radial 
modal spectrum does not have a structure similar to that of the highest peaks of 
the FS. Further, even if one insisted on nonradial modes, the visual observations are 
necessarily whole disk observations which cannot possibly resolve nonradial modes 
with £ > 2 . Again the required modal spectrum would not be compatible with stellar 
structure." 

The stars cannot be evolving (time-dependent) multiperiodic systems either. 
Otherwise there should be a correlation between the multitude of smaller peaks in 
the Fourier spectrum when subsequent sections of the data are Fourier analyzed, 
which is not the case. 

We conclude that neither the R Set nor the AC her lightcurve can reasonably 
be explained as being periodic or multi-periodic. 

3. The global flow reconstruction 

We now make the assumption that the observed lightcurve s n = {s(t n )} has 
been generated by a low dimensional dynamics in some 'physical' phase-space of a 
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priori unknown properties, i.e. 



dY/dt = G(Y) (1) 

where Y(t) e R d is the d dimensional state vector of the system. We could equivalently 
assume a stroboscopic description of the dynamics 

Y n+1 = G(Y n ) (2) 

where Y n = Y(t n ). It is natural to assume that s(t n ) is a smooth function of the 
phase-space variable Y. 

A standard reconstruction method starts with the construction of the standard 
<i e -dim delay vectorsiJl 

X n = { } (3) 

The global flow reconstruction method then assumes that there is a single function 
(map) F that describes the dynamics all over, 

X n+1 = F(X n ) (4) 

Embedding theorems! guarantee that for sufficiently large d e there is a diffeo- 
morphism between the physical phase-space and the reconstruction space (the latter 
is then called an embedding space, stricto sensu). Following Giona et aJ.i and BrownS 
we assume that the function F(X) is of a multivariate polynomial form, 

F(X) = 5>P. 4 (X) (5) 

i 

where the sum runs over all polynomials up to degree p. More specifically we use 
natural polynomials, i.e. polynomials that are constructed to be orthogonal on the 
data set. However, instead of constructing these polynomials explicitly by recursion 
and then determining the coefficients a, HI it is equivalent and much faster to expand 
the function in terms of all monomials up to the same order p 

F(X) = $>M(X) (6) 

i 

and then find the expansion coefficients fa through a least squares minimization of 
the variance 

s = J2 lix n - E AM(x n _i)|| 2 . (7) 

n i 

For this purpose a SVD (singular value decomposition) approach has been found to 
be most stable and efficient .LLiHJ 



The global polynomial reconstruction method was first tested^ on the Rossler 
attractor.El Using as input a short time-series of only the first one of the 3 Rossler 
variables one obtains a very good estimation of the Lyapunov exponents and the 
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fractal dimension of the attractor. Furthermore it is possible to determine that this 
scalar signal was generated by a 3-D flow (Rossler's three first order ODEs), in other 
words it is possible to determine not only the lowest embedding dimension, but also 
the dimension of the actual physical phase-space of the dynamics. 

Encouraged by these benchmark results we have applied the reconstruction 
method to the observational data of stellar lightcurves shown in Figs. 2 and 3. 

Because the observational data sets are so small we have used all the points 
in constructing the map in Eq. 7. We note in passing that in Astronomy we are not 
interested in prediction, but rather in uncovering the properties of the dynamics that 
underlies the pulsation. In Fig. 6 we display the behavior of the error as a function of 
d e , which therefore is an in-sample error. The vertically stacked points are for values of 
p ranging from 3 to 6. This suggests strongly that for R Set the minimum embedding 
dimension is d e = 4. For AC Her, because of the much larger observational noise, the 
lowest achievable error level is seen to be much higher. Nevertheless Fig. 6 can be 
seen to indicate that we may get away with a minimum d e of 3. The behavior of the 
error norm with dimension by itself is of course insufficient to establish the minimum 
embedding dimension. 
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Fig. 6 Behavior of the in-sample error for R Set and AC Her. 



It is worthwhile noting that the results obtained with the 'false nearest neighbor' 
methodi are consistent with these respective minimum dimensions of 4 and 3. 

How do we decide that a map is 'good' ? For that we have to compare the 'syn- 
thetic' signals that a map yields to the original data. (By synthetic signal we mean 
data that have been generated by iteration from some seed values with a map gener- 
ated according to Eq. 7). Of course, since we have chaotic behavior this comparison 
has to be done in an average or statistical sense. Unfortunately the observational data 
set is too small to make statistical comparison tests very meaningful. Instead we rely 
on three rather ad hoc criteria, viz. (a) the signals must have the same appearance, 
(b) their Fourier spectra must have the same envelope structure and (c) their lowest 
BK projections should be similar. We have found these simple criteria to be quite 
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useful. For example in our flow reconstructions it is quite clear that any 3D maps 
for R Set can be rejected for not being able to capture the dynamics of the star; the 
synthetic lightcurves, while they may be chaotic, they bear essentially no resemblance 
to the observational data. 

In Figs. 2 and 3, bottom rows, we exhibit the synthetic signals we could produce 
from what we consider the best of these maps for R Set and AC Her, respectively. In 
the BK projections of Fig. 6 we display the same two signals, Cols. 2 and 7, respec- 
tively. For comparison we also show synthetic signals for R Set that were obtained 
with a 5D map (Col. 3). Col. 6 shows a 3D reconstruction of AC Her. For a comparison 
of the Fourier spectra we refer to the original papers. 

For further details concerning the sensitivity to the delay A, the maximum 
polynomial order p and the data preparation can be found in the original papers.@! 

The reconstructions are seen to be very good and robust for R Set. Col. 4 of 
Fig. 6 shows that this remains true when we reconstruct the dynamics with a 4D flow 
instead of the 4D map of Col. 2. 

The AC Her data are inferior. Their signal to noise ratio is much lower and as a 
result the reconstruction is much less robust. It is nevertheless interesting that we are 
able to reconstruct a map that can produce synthetic signals with similar properties 
to the observational data. 

It is not possible to extract Lyapunov exponents and a fractal dimension directly 
form the observational data because of their short length. However we can obtain 
these quantities once we have constructed the map or flow. In Table 1 we show the 
Lyapunov exponents {A^}, ordered from largest to smallest, and of the Lyapunov 
dimension for some of the reconstructions for both R Set and AC Her. The latter 
is defined as 

K 

d L = K + J2WK + i\ (8) 

t=i 

where K is the largest value such that the sum is positive. The coefficients have been 
obtained from the maps that have been constructed with the listed parameters d e , A 
and p. 

4. Results 

We only note the most important points: 

(1) While there is some scatter in the exponents and in the associated fractal dimen- 
sion di, for all successful reconstructions di < 4 for R Set and 3 for AC Her, inde- 
pendently of the embedding dimension. This is an indication that we have achieved 
an embedding of the dynamics. 

(2) The largest of the Lyapunov exponents is always positive. This is a confirmation 
that the signal is indeed chaotic. 

(3) We have not shown the second Lyapunov exponent because it is always close 
to zero. This result is reassuring, because it indicates that we are actually close to 
describing a flow for which one of the exponents would be exactly zero (autonomous 
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Table 1. Lyapunov exponents and Lyapunov dimensions for R Scuti and AC Herculis 



R Scuti AC Herculis 
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0.0033 


-0.034 




2.10 


4 5 


4 
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0.0015 


-0.029 




2.05 



system) . 

(4) For R Set | A3 1 < Ai, whereas the inequality is reversed for AC Her, indicating a 
larger dissipation in the dynamics of AC Her. 

As physicists we would like to infer more than just Lyapunov exponents and 
fractal dimensions, because as such they are not very useful. Thus it is interesting 
that from the linearization of the R Set map about its fixed point we foundi two spiral 
roots, one unstable and the second stable with approximately twice the frequency. 
This led to the interesting physical interpretation that the chaotic dynamics is gen- 
erated by the nonlinear interaction of two radial modes of pulsation that are in a 
close resonance condition. Such a physical picture is supported by further theoretical 
considerations. As for AC Her, it is quite likely that two strongly interacting vi- 
brational modes also underlie this dynamics, but that the stronger dissipation alluded 
to just above compresses the dynamics down to 3D, in the same way the 2D Henon 
map shrinks to a ID logistic map for small values of b. 

5. Conclusions 

The flow reconstruction has been applied to real data, viz. the light-curve of the 
irregular variable stars R Scuti and AC Herculis. Not only are these data contaminated 
with a good amount of observational and other extrinsic noise, but we have no a priori 
information about the dynamics that has generated them. The analysis has shown 
that the dynamics which generates the light-curve of R Scuti has a dimension of four, 
i.e. the signal can be generated by a flow in 4-D. Furthermore the attractor has a 
fractal dimension « 3.1. For the star AC Her we find a phase-space dimension of 3 
and a fractal dimension 2.05^^^.2.5. 

There are no observational data sets available on the even less irregular W Vir 
stars. However, numerical hydrodynamical simulations and their analysis with the 
global flow reconstruction^! indicate also a 3-D phase-space and a fractal dimension 
of cLl ~ 2.01 - 2.05. A systematic trend seems to exist, viz. that as we go from the 
low luminosity, weakly irregular W Vir stars to the high luminosity, strongly irregular 
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RV Tau stars, the fractal dimension gradually increases. 

The fact that these stars have such a low dimensional dynamics must be consid- 
ered remarkable. They undergo large luminosity fluctuations, up to a factor of 40 for 
Set, with strong shock waves and ionization fronts criss-crossing the stellar envelope. 
Yet, the overall pulsational behavior is found to be relatively simple. 

As a result of this nonlinear analysis the nature of the irregular variability of 
these types of stars is no longer a mystery, but it has been shown to be simply the 
manifestation of an underlying low dimensional chaotic dynamics. 
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